Ustawienia

knitr::opts_chunk$set(echo = TRUE)
set.seed(23)

Biblioteki

library(dplyr)
library(ggplot2)
library(caret)
library(lattice)
library(reshape)
source("https://raw.githubusercontent.com/briatte/ggcorr/master/ggcorr.R")

Wczytanie danych

#wczytanie danych ze strony do data.frame
data_set <- read.csv("http://www.cs.put.poznan.pl/dbrzezinski/teaching/zed/elektrownie.csv")
#dodanie do data.frame kolumny zawierającej godzinę i miesiąc, wyodrębione z kolumny "data"
data_set$godzina <- as.numeric(format(as.POSIXct(factor(data_set$data),format="%m/%d/%Y %H:%M"),"%H"))
data_set$miesiac <- as.numeric(format(as.POSIXct(factor(data_set$data),format="%m/%d/%Y %H:%M"),"%m"))

Kod przetwarzający dane

not_empty_radiation <- data_set %>% filter(irradiamento != 0)
empty_radiation <- data_set %>% filter(irradiamento == 0)

for(i in 1:nrow(empty_radiation)) 
{
    row <- empty_radiation[i,]
    m <- not_empty_radiation %>% filter(miesiac == empty_radiation[i,"miesiac"], godzina == empty_radiation[i,"godzina"], idsito == empty_radiation[i, "idsito"])
    if (nrow(m) ==0)
    {
      empty_radiation[i,"irradiamento"] <- 0
    }
    else
    {
      empty_radiation[i,"irradiamento"] <- mean(m[["irradiamento"]]) #summarize(m, irra_mean = mean(irradiamento))
    }
} 
#połączenie dwóch data frameów, w wynikowy "data1"
data_set1 <- rbind(empty_radiation, not_empty_radiation)

Analiza danych

Zaimportowany zbiór zawiera dane pochodzące z czujników umieszczonych przy panelach fotowoltaicznych, umieszczonych w trzech sąsiadujących elektrowni słonecznych we Włoszech. Zbiór składa się z 235790 wierdzy oraz 51 kolumn. Celem analizy jest znalezienie czynnika, który najtrafniej pozwoli przewidzieć wartość wytworzonej energi przez panele fotowoltaiczne. Każdy wiersz zawiera uśrednione dane z jednej godziny pomiarów pojedynczej jednostki fotowoltaicznej.

Kolumny zbioru danych dotyczą specyfikacji panelu (identyfikator panelu, modelu, firmy), wartości czynników pogodowych (napromieniowanie, ciła wiatru, zachmurzenie, ciśnienie, wilgotność, punkt rosy), położenia geograficznego(długość i szerokość geograficzna panelu, wysokość, azymut), daty pomiaru (dzien, rok, pelna data) oraz wartości współczynnika “pcnm”.

W zbiorze danych znajują się wartości zerowe, nie tylko w miejscach gdzie żeczywiście odnotowany został taki pomiar, ale także gdzie, poprzez usterke panelu nie doszło do pomiaru. Powyższy wniosek opiera się na analizie pomiarów dokonanych o zbliżonej porze dnia, roku oraz stopniu zachmurzenia.

Położenie geograficzne nie jest bardzo istotne, ponieważ wszystkie panele są rozmieszczone w stosunkowo niewielkich odległościach (sasiednie elektrownie na terenie jednego kraju, takie samej strefy klimatycznej i wysokości nad poziomem morza). W kategorii czynników pogodowych, najbardziej istotne jest napromieniowanie, ponieważ to od stopnia nasłoneczenie zależy wartość kwh. Z tym bezpośrednio związana jest data pomiaru, ponieważ siła napromieniowania jest różna w zależności od pory dnia i roku. Dlatego też wartości “puste” najbardziej dotlikwe są dla danych dotyczących napromieniowania.

Nastepnie sporządzone zostały hostogramy dla atrybutów oraz korelacja pomiędzy poszczególnymi zmiennymi zbioru danych. Zgodnie z początkową analizą korelacja ukazała, ze na ilość wytwarzanej energii ma wpływ promieniowanie, a także wilgotność. Te argumenty zostały wskazane do regresji.

Sekcję podsumowującą rozmiar zbioru i podstawowe statystyki

#uwtorzenie podsumowania, statystyk ze zbioru danych "dataset". 
selected_data <- select(data_set1,idsito, idbrand, idmodel,temperatura_ambiente, irradiamento, pressure, windspeed, humidity, cloudcover,dewpoint, kwh, azimuth)
summary(selected_data)
##      idsito          idbrand          idmodel       temperatura_ambiente
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   Min.   :0.0450      
##  1st Qu.:0.1000   1st Qu.:0.0830   1st Qu.:0.1670   1st Qu.:0.2120      
##  Median :0.2250   Median :0.1670   Median :0.2080   Median :0.3480      
##  Mean   :0.2147   Mean   :0.1519   Mean   :0.2426   Mean   :0.3734      
##  3rd Qu.:0.3250   3rd Qu.:0.1670   3rd Qu.:0.2920   3rd Qu.:0.5300      
##  Max.   :0.4250   Max.   :0.4170   Max.   :0.7500   Max.   :0.8180      
##   irradiamento        pressure        windspeed          humidity     
##  Min.   :0.00000   Min.   :0.0000   Min.   :0.00000   Min.   :0.1600  
##  1st Qu.:0.00000   1st Qu.:0.7480   1st Qu.:0.04200   1st Qu.:0.5400  
##  Median :0.04273   Median :0.7530   Median :0.06600   Median :0.7000  
##  Mean   :0.11231   Mean   :0.6504   Mean   :0.07622   Mean   :0.6844  
##  3rd Qu.:0.20900   3rd Qu.:0.7550   3rd Qu.:0.10200   3rd Qu.:0.8400  
##  Max.   :0.71000   Max.   :0.7690   Max.   :0.69600   Max.   :1.0000  
##    cloudcover       dewpoint           kwh            azimuth      
##  Min.   :0.000   Min.   :0.1390   Min.   :0.0000   Min.   :0.1280  
##  1st Qu.:0.230   1st Qu.:0.5350   1st Qu.:0.0000   1st Qu.:0.2950  
##  Median :0.310   Median :0.6190   Median :0.0490   Median :0.4250  
##  Mean   :0.359   Mean   :0.6055   Mean   :0.1688   Mean   :0.4546  
##  3rd Qu.:0.510   3rd Qu.:0.6830   3rd Qu.:0.3320   3rd Qu.:0.6350  
##  Max.   :1.000   Max.   :0.8650   Max.   :1.0000   Max.   :0.8180

Szczegółowa analizę wartości atrybutów

#Wybranie określonych kolumn ze zbioru dataset, po wcześniejszej pobieżnej analizie ich ważności. Wybrane zostały głównie kolumny związane z czynnikami pogodowymi a także identyfikator panela fotowoltaicznego oraz wartość prezentującą wytworzoną energię. 
selected_data <- select(data_set1,idsito, idbrand, idmodel,temperatura_ambiente, irradiamento, pressure, windspeed, humidity, cloudcover,dewpoint, kwh, azimuth)

#W poniższym kroku konstruowany jest wykres za pomocą funkcji "ggplot", który przedtawia rozkład wartości każdego z aktrybutów. 
p <- ggplot(data = melt(selected_data), aes(x = value)) + 
  geom_histogram(bins=10, aes(fill=..count..)) + 
  labs(title = "Rozkład wartości atrybutów - histogram") + 
  facet_wrap(~variable, ncol=4) + 
  scale_x_continuous(labels = scales::comma) + 
  scale_fill_gradient("Count", low="blue", high="pink") +
  theme_bw()
## Using  as id variables
ggplotly(p)

Korelacja pomiędzy wybranymi atrybutami

selected_data <- select(data_set1,idsito, idbrand, idmodel, ageinmonths, temperatura_ambiente, irradiamento, pressure, windspeed, humidity, cloudcover,dewpoint,godzina, miesiac, kwh)
ggcorr(data_set, nbreaks = 4, label = TRUE, label_size = 2, label_color = "white")
## 
## Attaching package: 'reshape2'
## The following objects are masked from 'package:reshape':
## 
##     colsplit, melt, recast
## Warning in ggcorr(data_set, nbreaks = 4, label = TRUE, label_size = 2,
## label_color = "white"): data in column(s) 'data' are not numeric and were
## ignored

Prezentacja zmiany wytwarzanej energii w czasie i przestrzeni

data_set <- data_set %>% mutate(data_format=format(as.POSIXct(data, format='%m/%d/%Y %H:%M'), "%Y-%m"))
kwh_plot <- data_set %>% group_by(data_format, idsito) %>% summarise(suma=sum(kwh))
p <- ggplot(data = kwh_plot, aes(x=data_format, y=suma, color=factor(idsito))) + 
  geom_point() + labs(title = "Wykres energii w czasie i przestrzeni", x = "Data", y = "Ilość wytworzonej energii[kwh]") + theme(axis.text.x = element_text(angle = 90))

ggplotly(p)

Regresor przewidujący wytwarzaną energię

set.seed(23)
df_clear <- data_set %>% select(idsito, irradiamento, humidity, kwh)

inTraining <- createDataPartition(
    y = df_clear$idsito, 
    p = 0.85,
    list = FALSE
)

training <- df_clear[inTraining,]
testing <- df_clear[-inTraining,]

ctrl <- trainControl(method="repeatedcv", number=2, repeats=5)

model <- train(
    kwh ~ idsito + humidity + irradiamento,
    data = training,
    method = "lm",
    metric = "RMSE",
    trControl = ctrl
)

predict_kwh <- predict(model, testing) 
dat <- data.frame(pred =  predict_kwh, obs = testing$kwh)
defaultSummary(dat)
##       RMSE   Rsquared        MAE 
## 0.10639435 0.74410469 0.06410652

Analizę ważności atrybutów

Na podstawie analizy atrybutów oraz modelu regresji można stwierdzić, że do przewidzenia wartość wytworzonej energiiprzez pojedynczy panel w danej godzinie najlepiej posłuży parametr dotyczący napromieniowania, wilgotnosci oraz azymutu. Trafność regresji została oszacowana na podstawie miary RMSE i jej wartość jest na poziomie 0.10.

#Podsumowanie modelu
print(summary(model))
## 
## Call:
## lm(formula = .outcome ~ ., data = dat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.66577 -0.03691 -0.01354  0.01887  1.01985 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   0.075492   0.001279   59.03   <2e-16 ***
## idsito        0.095897   0.001776   53.99   <2e-16 ***
## humidity     -0.098289   0.001558  -63.09   <2e-16 ***
## irradiamento  1.282707   0.002109  608.29   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1056 on 200418 degrees of freedom
## Multiple R-squared:  0.7487, Adjusted R-squared:  0.7487 
## F-statistic: 1.991e+05 on 3 and 200418 DF,  p-value: < 2.2e-16